Métodos Funcionales en Aprendizaje Automático

Examen 02 - Manifold Learning

José Antonio Álvarez Ocete

Configuración Inicial

Directrices generales

Este notebook contiene el enunciado de la segunda parte del examen de Métodos Funcionales en Aprendizaje Automático.

Por favor, lee con cuidado cada uno de los tres ejercicios así como las preguntas que en ellos se formulan, y explica de manera concisa cada asunción y conclusión que hagas.

Debéis entregar un informe (que puede realizarse directamente sobre el notebook) y el código implementado para completar los ejercicios, así como las referencias consultadas en caso de necesidad. También se debe entregar el código de honor firmado. Todo esto se enviará vía Moodle en un .zip.

Ejercicio 1 (3 puntos)

grafo.png

* Para el grafo anterior, calcula $f^\top L f$, siendo $L$ el laplaciano del grafo no normalizado ($L = D- W$) cuando $f$ se define como: $f(x_1) = 1, f(x_2) = 2, f(x_3) = 6, f(x_4) = 5, f(x_5) = 4$ si los pesos del grafo son: $w_1 = 20, w_2 = 10, w_3 = 1, w_4 = 5$.

Solución.-

Asumimos que el peso de los lazos (lados conectado un nodo consigo mismo) es 0, al no aparecer ningún peso en la imagen anterior. Definimos la matriz $W$ y computamos el laplaciano del grafo:

* Calcula $f^\top L f$ si cambiamos el valor de los pesos $w_1=1$ y $w_3=20$. ¿Cambia el resultado? ¿Por qué?

Solución.-

Como podemos ver, el resultado ha cambiado de 210 a 666. Para responder a por qué ocurre este hecho hemos de dar contexto a las operaciones que estamos realizando. En el contexto de los Laplacian Eigenmaps, se crea un grafo para un conjunto de datos con pesos asociados a la afinidad de los nodos. Es decir, cuanto más cerca estén dos muestras $x_i, x_j$, mayor será el peso asociado entre las mismas, $w_{ij}$.

En las notas de clase definimos el error de reconstrucción:

$$ \mathcal J(y) = \frac{1}{2} \sum_{i,j}(y_i - y_j)^2 w_{ij}. $$

Para minimizar este error hemos de asignar valores de $y_i$ parecidos a elementos con altos pesos (es decir, si los puntos son cercanos, asignamos valores de la función parecidos para que la reconstrucción tenga sentido). También demostramos en clase que el error de reconstrucción coincide con la evaluación de la siguiente forma cuadrática:

$$ \mathcal J(y) = y^T L y. $$

En nuestro caso particular, los valores de $f(x_1)$ y $f(x_2)$ son parecidos ($1$ y $2$ respectivamente), mientras que el valor de $f(x_3) = 6$ dista más de ellos. Por lo tanto, para minimizar el error de reconstrucción, los pesos han de cumplir $w_{1,2} \gg w_{2,3} > w_{1,3}$.

Inicialmente: $w_{1,2} = 20, w_{1,3} = 1$, lo que tiene sentido (los nodos 1 y 2 estarán muy cerca mientras que los nodos 1 y 3 estarán muy lejos). En el segundo caso $w_{1,2} = 1, w_{1,3} = 20$, representando a los nodos 1 y 2 muy lejanos y los nodos 1 y 3 muy cercanos, lo que no encaja para nada con los valores de $f(x)$. Es por ello que el valor de $f^T L f$ (el error de reconstrucción) aumenta considerablemente.

* Define una función $f$ que minimice $f^\top L f$. ¿Existe más de una solución que minimice $f^\top L f$ , sin tener en cuenta múltiplos de la $f$ definida (de manera más formal, que sea linealmente independiente a la solución dada)?

Solución.-

Nota inicial: Realizaremos la minimización para ambas matrices $W$ definidas en los apartados anteriores. La solución obtenida es válida considerando cualquiera de las dos.

Sabemos que el laplaciano del grafo no normalizado ($L = D - W$) es una matriz semidefinida positiva. Esto es, $x^T L x \ge 0$ para todo vector $x$. El vector $f_1 = 0$ obtiene como resultado $f_1^T L f_1 = 0$, minimizando la función.

Puesto que nuestra función $g(x) = x^T L x$ es una forma cuadrática, sabemos que es analítica (en particular dos veces derivable), y que una función dos veces derivable es convexa si y sólo si su matriz Hessiana ($L$, por ser una forma cuadrática) es semidefinida positiva, y es estrictamente convexa si y sólo si $L$ es definida positiva.

Si $f$ fuese estrictamente convexa, tendríamos un único mínimo. Sin embargo, sabemos que $\lambda_0 = 0$ siempre es valor propio de $L$, como vimos en clase. En nuestro caso particular:

Al tener un valor propio $0$ sabemos que $L$ es semidefinida positiva, luego $g$ no es estrictamente convexa y puede tener otro valor mínimo. Para buscarlo intentamos en primer lugar vectores con valores no nulos en una única posición, pero esto nos devuelve valores de la diagonal multiplicados por un número positivo, por lo que nunca obtendremos un cero.

Hemos de hacer uso de los valores negativos, no situados en la diagonal. Mirando la submatriz inferior derecha de $L$ vemos que si utilizamos el vector $f_2 = (0, 0, 0, 1, 1)$ obtendremos un valor nulo, ya las dos últimas filas (y por tanto las dos últimas columnas) suman $0$ cada una:

Sin embargo, con esto no basta para asegurar que existan soluciones linealmente independientes, pues $f_1$ y $f_2$ son linealmente dependientes: $f_1 = 0 \cdot f_2$. Mirando la submatriz superior izquierda del laplaciano sin normalizar con la idea anterior en mente vemos que la suma de las tres primeras filas (así como la de las tres primeras columnas) es nula, por lo que también nos sirve la solución $f_3 = (1, 1, 1, 0, 0)$:

Con esto hemos encontrado dos mínimos de $g$ linealmente independientes: $f_2$ y $f_3$.

Ejercicio 2 (4 puntos)

Un algoritmo muy utilizado tanto en el ámbito de manifold learning, como en clustering o en análisis de grafos (por ejemplo para detectar comunidades) es el método de Spectral Clustering. Durante las clases del curso no nos ha dado tiempo a profundizar en él, por lo que vamos a utilizarlo para resolver un problema.

* Define una clase **Spectral Clustering**, así como sus funciones fit y fit-transform (no vamos a permitir en este caso la aplicación del método a puntos nuevos fuera de la muestra inicial). * El algoritmo de Spectral Clustering a implementar será el basado en el laplaciano de camino aleatorio, es decir, $L_{rw}=I − D^{−1}W$.

Solucion.-

El algoritmo de Spectral Clustering está compuesto por dos pasos clave:

La siguiente celda implementa dicho algoritmo de la forma más genérica posible y utilizando como base la clase de DiffusionMap creada para la práctica. Esta incluye:

  1. La matriz de pesos del grafo asociado, W, (también denominada matriz de afinidad) puede calcularse de distintas formas, elegida por el parámetro affinity:

    • (Por defecto) Si affinity == 'kernel', se utilizará un kernel para este computo. La posición $w_{ij}$ de la matriz será kernel(x_i, x_j). El kernel utilizado se define a partir del parámetro kernel, pudiendo ser rbf (por defecto), laplacian o una función kernel pasada por parámetro. La longitud del kernel ($\sigma$ o $\gamma$ dependiendo del kernel utilizado) se define a partir del parámetro sigma. Este sigma puede calcularse automáticamente utilizando los valores de median o maximum, como se explico en la práctica de DiffussionMap.

    • Si affinity == knn, se utilizará el algoritmo de $k$ vecinos más cercanos (KNN) para computar la matriz de afinidad. En particular $w_{ij}$ será $1$ si el elemento $x_j$ está entre los $k$ vecinos más cercano del elemento $x_i$. Como está relación no tiene por qué ser simétrica, tras el cálculo de la $W$ inicial la simetrizamos haciendo la media con su traspuesta. El número de vecinos queda definido por el parámetro n_neighbors.

    • Si affinity == precomputed, la matriz de afinidad se ha calculado previamente y se suministra a través del parámetro precomputed_W.

  2. Una vez la matriz de pesos ha sido calculada, el método de obtención del Laplaciano (fijado con el parámetro laplacian_alg) varía entre:

    • Si laplacian_alg == unnormalized: utilizamos el laplaciano sin normalizar: $L = D - W$.
    • (Por defecto y el requerido para este examen) Si laplacian_alg == random walk: utilizamos el laplaciano normalizado por la izquierda, también denominado random walk: $L = I - D^{-1} W$.
    • Si laplacian_alg == normalized: utilizamos el laplaciano normalizado: $L = I - D^{-\frac{1}{2}} W D^{-\frac{1}{2}}$.

  3. Una vez hemos computado el laplaciano, lo diagonalizamos y proyectamos a un espacio con n_components dimensiones. Para ello nos quedaremos con los n_components vectores propios de $L$ con menor valor propio asociado. Si drop_first == True, no consideraremos el primero vector propio, por estar asociado al valor propio $0$.

    • Nota: en este punto he estudiado definir un threshold y eliminar todos los vectores asociados a valores propios con valor menor que dicho threshold. Sin embargo, en el dataset de S-curve hay varias decenas de valores propios con valor menor a $10^{-15}$, y adicionalmente Sklearn no implementa esta disquisición, por lo que decidí no añadirla.

  4. Una vez obtenido el embedding, aplicamos el algoritmo K-means de Sklearn de forma directa, donde el número de clusters queda definido por el parámetro n_clusters.

Como nota final, todos los cálculos potencialmente aleatorios tienen la semilla fijada. Esta puede controlarse con el parámetro random_state.

Finalmente, para diseñar e implementar esta clase se han consultado las siguientes referencias:

[1]: Wikipedia - Spectral clustering.

[2]: Sklearn - Spectral clustering.

[3]: Github - Sklearn's spectral clustering implementation.

[4]: Sklearn - Nearest Neighbors.

* ¿Qué diferencias aprecias entre este algoritmo y *Diffusion Maps*?

En primer lugar, Spectral Clustering es un algoritmo de clustering, no de embedding. Podemos comparar Diffusion Maps con el algoritmo que utiliza Spectral Clustering para obtener el embedding. Esto es, Laplacian Eigenmaps con distintos métodos para crear el laplaciano.

Por un lado, ambos métodos buscan utilizar la geometría intríseca de los datos para reducir la dimensionalidad de los mismos, de forma no supervisada. Sin embargo, la motivacición de los métodos es originalmente completamente distinta: Diffusion maps busca simular un proceso de expansión del calor para estudiar dicha geometría, mientras que Laplacian Eigemaps construye un grafo y busca minimizar el error de reconstrucción, por lo que transforma el problema en minimizar un Lagraniano, aunque llegue a una solución similar en última instancia.

Respecto a cómo llegar a dicha solución, los algoritmos difieren entre si en los siguientes detalles:

Respecto a los parámetros de cada algoritmo, ambos necesitan la dimensionalidad a la que proyectar, así como los detalles de cómo computar la matriz de afinidad (el kernel utilizado si es que se utiliza un kernel, así como los parámetros del mismo). Sin embargo, Difussion Maps tiene un parámetro adicional, el número de pasos a simular. Esto lo hace algo más complicado de ajustar.

* Aplica el algoritmo definido al siguiente conjunto de datos. * ¿Se obtiene un buen resultado? * ¿Cómo has realizado la selección de hiperparámetros involucrados?

Para responder a las preguntas vamos a hacer algunos experimentos con el dataset proporcionado.

Embeddings

Cargamos el dataset proporcionado

En primer lugar mostramos únicamente el embedding obtenido utilizando Laplacian Eigenmaps proyectado a dos dimensiones, asi como las dos primeras dimensiones (de doce) de los datos originales. Utilizamos el kernel rbf con sigma la mediana para obtener el embedding.

Como podemos apreciar, los datos se separan en algo parecido a dos medias lunas intersecadas. Como hemos estudiado amplicamente, este dataset es relativamente sencillo de separar utilizando algún método no lineal no tirivial (como una SVM con kernel no lineal). Sin embargo, es obvio que a la hora de hacer clustering no obtendremos clusters que encaja particularmente bien con las clases originales.

Volvemos a ejecutar este mismo ejemplo utilizando el algoritmo de KNN para computar la matriz de afinidad con 10 vecinos (el valor por defecto).

En esta ocasión los puntos obtenidos son algo más separables linealmente, aunque la intersección no es trivial y hay varios puntos azules dentro del clúster rojo. Intuitivamente este método funcionará mejor a la hora de hacer clustering.

Motivados por la segunda pregunta del enunciado, realizamos un 'mini-gridsearch' sobre valores de $\sigma$ y $k$ para los métodos anteriores, elegidos a mano. Mostramos a continuación los resultados:

Para el primer experimento sobre posibles valores de $\sigma$, vemos como a partir de cierto valor de sigma ($\sim 1$) el embedding se estabiliza y ya no varía. Las dos útlimas gráficas son las del la mediana y el máximo respectivamente. Como se puede apreciar, el embedding obtenido por los valores por defecto es esencialmente el mejor.

Respecto al segundo experimento, vemos como el embedding no llega a estabilizarse para los valores de $k$ probados. Sin embargo, conforme se aumenta lo suficiente, los clusters empiezan a mezclarse significativamente. Se puede apreciar en los valores de $k$ mayores de $10$. Sin embargo, el caso de $k=50$ separa relativamente los clusters, aunque tampoco significativamente mejor que en el caso de RBF.

Como sabemos, los algoritmos de manifold learning tienen como objetivo aprovechar la geometría original de las variedades diferenciales en las que se encuentran los datos para desenrollarlas en menor dimensión. Por lo que estamos observando, no parece que los datos estén contenidos en una variedad diferencial desenrollable en dos dimensiones. Probemos con tres:

Tampoco parecen proyectables en una variedad diferenciable tres-dimensional.

Clustering

A continuación probamos el algoritmo completo de Spectral Clustering, que incluye el paso de obtención de clusters tras computar el embedding.

Como ya predijimos, el clustering no separa las clases en clusters, ya que estaos no separan particularmente bien los datos. Si quisiesemos clasificar los datos deberíamos aplicar alguna técnica de aprendizaje supervisado no lineal, que debería poder separar las clases de forma mucho más sencilla.

Preguntas planteadas

Respondemos a las preguntas enunciadas al principio de este ejercicio.

El embedding nos ayuda a separar las clases mejor que si no lo utilizásemos, pero desde luego no separa los datos de forma obvia.

En primer lugar he probado con los parámetros por defecto, tanto los de Sklearn ($\sigma = 1 / X.shape[1]$) como los computados en la práctica anterior (el máximo y la mediana). De la misma forma, para el método de KNN comenzamos utilizando el valor por defecto de Sklearn, $k=10$.

Viendo que estos parámetros no separan claramente las clases, he realizado un grid-search manual para estudiar otros posibles valores de $\sigma$ y $k$, sin encontrar mejores resultados significativos.

S-curve dataset

Comparación con Sklearn

Para garantizar que mi implementación de Spectral Clustering es correcta, utilizaremos el conocido dataset de la curva en S tridimensional e intentaremos desenrollarlo en dos dimensiones. Comenzamos cargando los datos:

A continuación utilizaremos el algoritmo de Spectral Clustering implementado por Sklearn y el nuestro sobre los mismos datos para ver si obtenemos el mismo resultado, tanto para el embedding obtenido como los clusters.

En primer lugar cabe destacar que la expresión del kernel RBF y la nuestra difieren en la colocación del parámetro:

$$ f_{sklearn}(x,y) = e^{- \gamma \| x - y\|^2} \quad \quad f_{ours}(x,y) = e^{- \frac{1}{2\sigma^2} \| x - y\|^2} $$

Utilizando el valor de parámetro por defecto de Sklearn ($\gamma = 1$) hemos de utilizar $\sigma = \frac{1}{\sqrt 2}$.

Vemos como los embeddings obtenidos son muy parecidos pero no son ambos buenos y realmente parecidos, pero no encajan exactamente. Comentaremos este hecho más adelante.

Respecto a los clusters obtenidos, ambos algoritmos obtienen resultados realmente buenos: los clusters corresponden a secciones de la superficie desenrollada, y no a clusters obtenidos utilizando una distancia euclídea en el espacio original tridimensional. Este era el resultado esperado por nuestro algoritmo y confirma que la implementación es correcta.

En el caso del clustering realizado por Sklearn vemos como en los límites de la S, los puntos se asignan a un cluster distinto, no al más cercano en la superficie. Esto cobra sentido mirando el embedding obtenido, donde los límites de la S quedan relativamente cercanos a la zona central. Si colocamos un centro de cluster en esta zona central, los puntos del límite se asignarán a este clúster, como ocurre en nuestro caso.

Finalmente, comentemos por qué no obtenemos exactamente el mismo embedding (y por tanto el mismo clustering) respecto al obtenido por Sklearn. Tras mucho estudio sobre este hecho, la conclusión final es que Sklearn no implementa exactamente el mismo algoritmo de clustering que nosotros. Para aclarar este hecho he tenido que bucear en la implementación de Sklearn:

Adicionalmente, cabe destacar que para obtener el embedding realizado por el algoritmo de spectral clustering de Sklearn hemos de obtener a mano su matriz de afinidad y los vectores propios. Sin embargo, no hemos de utilizar los asociados a los menores valores propios sino aquellos asociados a los mayores valores propios, lo que ratifica que la implementación del embedding es distinta.

Embeddings

A continuación realizamos un experimento similar al realizado para el dataset proporcionado, mostrando distintos embeddings para distintos valores de los parámetros, tanto para $\sigma$ como para $k$, para observar si podemos mejorar el embedding obtenido.

Por un lado, utilizando el kernel RBF para computar la matriz de afinidad obtenemos formas de $V$ en nuestro embedding que despliegan relativamente bien nuestra superficie (sin desplegarlo en forma de alfombra como un rectángulo perfecto, que sería lo óptimo intuitivamente). A partir de cierto valor de $\sigma$ ($\sim 1$) obtenemos una forma de $S$, lo que captura muy bien la informaciónde de la superficie enfocada desde un lateral, pero no la extiende bien, que sería nuestro objetivo.

Por otro lado, utilizando KNN para computar la matriz de afinidad obtenemos un resultado relativamente parecido mientras aumentamos el número de vecinos utilizados. Sin embargo, las formas de $V$ obtenidas parecen extender peor la superficie a lo ancho para valores bajos de $k$. Además, si bien nos acercamos a obtener esa forma de $S$ como ocurre en RBF, esta es mucho menos consistente y la forma queda menos definida.

Ejercicio 3 (3 puntos)

Hasta ahora hemos visto cómo trabajar con una matriz de datos que transformamos en un grafo para aplicar un método de reducción de dimensión espectral, pero todos los algoritmos que hemos vistos podrían aplicarse directamente a un grafo con cierta información. Para este ejercicio vamos a trabajar con el grafo del club de karate de Zacarías (un ejemplo clásico en teoría de grafos).

Para leer los datos necesitaréis instalar la librería igraph y cairo, para lo que podéis usar los comandos: pip install python-igraph pip install pycairo Y si os da error porque faltan librerías deberéis instalar: sudo apt install libcairo2-dev sudo apt install python3-dev

Veamos primero que pinta tienen los datos:

El siguiente código permite visualizar el grafo dado. En este caso permitimos que se seleccione de manera automática el mejor algoritmo de visualización, etiquetamos los nodos por su identificador y lo coloreamos (todo del mismo color pues a priori no conocemos las comunidades).

* Obten la matriz de pesos del grafo a partir de los datos dados usando la librería igraph proporcionada.

Puesto que los lados de este grafo no tienen pesos asociados (la forma más sencilla de apreciar esto es observar los datos en el documento proporcionado, karate.gml), la matriz de pesos del grafo será la matriz de adyacencia del mismo. Para obtenerla, recorremos los lados del grafo y y añadimos las aristas correspondientes a nuestra matriz.

* Aplica el algoritmo de Spectral Clustering de manera que se observen lo mejor posible las dos comunidades existentes en el grafo.

Creamos una función que divide el grafo elegido (a partir de la matriz de pesos pasada por parámetro, W), en tantos clusters como elijamos (n_clusters) y proyectando a tantas dimensiones como queramos (n_components). Variando el parámetro method podemos elegir si utilizamos la implementación de Spectral Clustering de Sklearn o la nuestra.

Mostramos en primer lugar la nuestra y en segundo lugar, la de Sklearn.

Vemos como nuestra aproximación (reduciendo a dos dimensiones y tomando dos clusters) no captura tan bien como Sklearn la comunidad inferior derecha. Explorando con distinto número de dimensiones he encontrado que la mejor dimensión a la que proyectar para obtener ambas comunidades es n_components=1, siendo casi equivalente a a la de Sklearn pues sólo se diferencian en un nodo, el $20$.

Aunque puede parecer contraintuitivo proyectar en una única dimensión antes de realizar el clustering, esto es esencialmente proyectar en una recta y separar los dos clusters cortando la recta en algún punto.